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ABSTRACT 

We present continuum observations at A = 1.3 mm and 2.7 mm using the 
Combined Array for Research in Millimeter-wave Astronomy (CARMA) toward 
six protoplanetary disks in the Taurus molecular cloud: Cl Tan, DL Tan, DO 
Tan, FT Tan, Haro 6-13, and HL Tan. We constrain physical properties of the 
disks with Bayesian inference using two disk models; flared power-law disk model 
and flared accretion disk model. Comparing the physical properties, we hnd that 
the more extended disks are less flared and that the dust opacity spectral index 
(/3) is smaller in the less massive disks. In addition, disks with a steeper mid¬ 
plane density gradient have a smaller P, which suggests that grains grow and 
radially move. Furthermore, we compare the two disk models quantitatively and 
find that the accretion disk model provides a better £t overall. We also discuss 
the possibilities of substructures on three extended protoplanetary disks. 

Subject headings: circumstellar matter — protoplanetary disks — radio contin¬ 
uum; stars — stars: pre-main-sequence — techniques: interferometric 
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Introduction 


Circumstellar disks of young stellar objects (YSOs), particularly at the stages of Class 
II and III (the so-called T Tauri stars), are often called protoplanetary disks, since they are 
expected to form planets. These disk structures have been studied by spectral energy dis¬ 


tribution (SED) over the past few decades (e.g., Beckwith et ah 1990 Andrews & Williams 
2005). However, in order to derive the detailed physical properties radio interferometry 


at (sub)millimeter wavelengths is required, since the disks are cold (T ~ 30 K) and small 
{D < 1"). To date about two dozen protoplanetary disks have been studied by interfer¬ 
ometric observations with sub-arcsecond resolution (e.g, Andrews et ah 2009; Isella et ah 


2009 Guilloteau et ah 2011). Recently the Atacama Large Millimeter/sub millimeter Array 


(ALMA) started to present excellent results on transition disks, which are more evolved 
than protoplanetary disks, with the unprecedented sensitivity (e.g., van der Marel et ah 


2013 Carpenter et ah 2014). Also, detailed studies of grain growth in protoplanetary disks 


beyond millimeter sizes has been carried out with the expanded Karl G. Jansky Very Large 


Array (JVLA) (e.g., Perez et ah 2012). In addition, the long baseline science verification 
(SV) data of ALMA have recently shown unprecedented details of substructures toward HL 
Tan, which is one of our sample. 

In the previous modeling studies of protoplanetary disks using SED or early interfer¬ 
ometric data, a power-law disk model with sharp inner and outer edges was utilized (e.g., 
Beckwith et al.||1990 Mundy et al.||1996). Although the assumption is quite reasonable since 


the minimal-mass solar nebula shows a power-law density distribution (Weidenschilling 1977) 


and numerical simulations of circumstellar disks also show a power-law density distribution 


e.g., Ayliffe & Bate 2009), it has no fundamental physics behind it. In addition, the model 


does not explain the difference between disk sizes detected in dust continuum (compact) and 
GO spectral line (extended). 

The viscous accretion disk model provides a more physically motivated model, based 
on accretion with conservation of angular momentum (e.g., Pringle 1981). The model has 


a density distribution with a form of a power-law tapered by an exponential function. One 
complication in this model is that the sources of viscosity in the disk are not well understood 
and may vary with radius and disk mass. Hughes et al. (2008) argued that the viscous 


accretion disk model explains both dust continuum and gas spectral line data better than 
the power-law disk model by qualitative comparison. In this study, we employ both flared 
power-law disk and flared viscous accretion disk models to investigate disk properties; we 
compare the two disk models in order to determine which model better fits our millimeter 
data. Although some studies have modeled the dust grain properties with radius in the disk 
using submillimeter to centimeter data (e.g., Perez et al. 2012), we assume constant dust 
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properties in the disk due to our relatively limited wavelength coverage. 


Bayesian inference is used to determine the model parameter distributions that £t the 
observational data, instead of the commonly used method. The Bayesian approach pro¬ 
vides a basis of probability theory to compare two models as well as to achieve parameter 
probability distributions, which are not given by the method. Although the x^ method 
has been used widely, a couple of circumstellar disk studies have been carried out utilizing 
Bayesian inference ( Lay et ^|1997 Isella et ah [2009 ; Kwon||2009 ). No model comparison has 
been attempted for circumstellar disks so far in the Bayesian approach. Isella et ah (2010) 


obtained slightly smaller reduced x^ values toward their A = 1.3 mm data of DG Tan and 
RY Tan for the accretion disk model than the power-law disk model. On the other hand. 


Guilloteau et al. (2011) and Kitamura et ah (2002) also utilized both models. However, their 
studies used thin disk models and a simple power-law temperature distribution. Our disk 
models are flared and our temperature distributions are obtained by a Monte Garlo radiative 
transfer code, which are more realistic. 


In this paper, we present observational and modeling results of our protoplanetary disk 
survey, which has taken high angular resolution data at A = 1.3 mm and 2.7 mm using 
the Gombined Array for Research in Millimeter-wave Astronomy (GARMA) providing high 
image fidelity. First, brief introductions of the targets are given in Section and GARMA 
observations and data reduction are discussed in the following section. Afterward the two 
disk models of flared power-law disks and viscous accretion disks are explained, followed 
by observational and fitting results and discussion in Section Finally, we summarize the 
results. 


2. Targets 


The circumstellar disk targets of this study (GI Tan, DL Tan, DO Tan, FT Tan, Haro 
6-13, and HL Tan; Table are located in the Taurus molecular cloud, a well-known, nearby 


star forming region. The distance has been determined by various methods (Rebull et al. 


2004 

): e.g., star counting (e.g.. 

McGuskey||1938 

), optical extinction (e.g.. 

Kenyon et al.||1994 


140 ± 10 pc), parallax measurement (e.g., Hipparcos astrometric data: Bertout et al. 1999), 


and protostellar rotational properties (Preibisch & Smith 1997, 152 ± 10 pc). The measured 
distances to the Taurus molecular cloud are somewhat different in methods and in regions. 


For example, Bertout et al. (|1999|) reported three regions of 125l^g pc, IdO]*]}^ pc, and lOS^'^^s 


1+16 


3+42 


pc, using Hipparcos astrometric data. However, 140 pc is within uncertainties of most of 
the ranges, so it is adopted for all 6 targets in this study. Note that the distance affects the 
physical sizes and mass estimates. 
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Our targets were chosen from the sample of Beckwith et ah (1990), which observed 


86 YSOs in the Taurus molecular cloud at A = 1.3 mm using the IRAM 30 m telescope. 
They were also included in the submillimeter survey at A = 450 and 850 /rm carried out 


by Andrews & Williams (2005) using SCUBA on the James Clerk Maxwell Telescope. On 


the other hand, Kenyon & Hartmann (1995) compiled published IR and optical observations 


with supplementary IR observations in the cloud. In addition, Furlan et ah (2009) reported 
A ~ 10 /im feature strengths (equivalent width of the 10 /rm feature, EW(10 /rm)) and 
spectral index (ni 3 _ 3 i) between A = 13 and 31 /rm of T Tauri disks in the Taurus molecular 
cloud including our targets except HL Tau, as well as two other star forming regions using 
the Spitzer Space Telescope IRS data. They utilized these values to indicate dust evolution; 
for example, smaller ni 3 _ 3 i values represent disks with more dust settlement (less flaring) 
and EW(10 /im) shows the amount of small dust grains (< 5 /rm), which means that smaller 
EW(10 /im) indicates more evolved dust. Selected information from the above studies is 
summarized in Table More detailed studies have been carried out toward these targets 
and some selected examples are introduced in the following. 


Cl Tau was observed by [Beckwith et al.| ( |1990[ ) with a continuum flux of 190 ± 17 mJy 
at A = 1.3 mm. In addition, variability in optical spectra has been detected, which can be 


interpreted as a temporary obscuration of a local hot region (Smith et ah 1999). Andrews 


& Williams (2007) observed Cl Tau using the Submillimeter Array (SMA) at A = 880 p,m 
with a moderate angular resolution (> 1"). By model htting of SED data and visibility data 
averaged in annulus, they constrained physical properties such as density and temperature 


distributions, size, and mass. Guilloteau et ah (2011) also studied this target by about 
0.5" X 0.3" angular resolution data at two wavelengths of the Plateau de Bure Interferometer 
(PdBI) and two disk models: power-law and accretion disk models. 

DL Tau has been detected with a continuum flux of 230±14 mJy by IRAM observations 


(Beckwith et al.|1990 

. It has also been detected 

in the GO J = 2 —)■ 1 line, and the velocity 

held was shown 

(e-g-, 

Koerner & Sargent 

1995 

)• 

Simon et al. 

(2000 

) studied the kinematics 


using PdBI data in CO J = 2 —)■ 1. They showed that the disk is in Keplerian rotation and 
estimated the protostellar mass as 0.72±0.11 Mq assuming the inclination determined by the 
continuum image (49 ± 3°). Kitamura et ah (2002) has studied the disk structure by htting 
SEDs and Nobeyama Millimeter Array data taken at A = 2 mm with 1" resolution. They 


employed a viscous accretion disk as well as a power-law disk model. Andrews fc Williams 


(2007) also studied the DL Tau disk structure using SMA data at A = 1.3 mm with about 


2" resolution. This target was also studied by Guilloteau et al. (2011) using PdBI. 


DO Tau has a continuum hux of 136±11 mJy in the Beckwith et ah (1990) observations. 
Koerner & Sargent (1995) have reported GO J = 2 —)■ 1 observational results . The line 
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profile has a much broader feature in the blueshifted region than in the redshifted region, 
which indicates that it needs infall and bipolar outflow components. DO Tan is also one of 


the sample in the studies of Kitamura et ah (2002). 


FT Tau has a continuum flux measured in Beckwith et ah (1990) of 130 ± 14 mJy at 


A = 1.3 mm. In previous survey observations at optical wavelengths, however, the spectral 
type of FT Tau was not determined probably due to an undetermined luminosity, although 
the surveys included the target. Therefore, there is no good estimate of protostellar mass 
and age, which are normally determined by its effective temperature and luminosity. Ro- 


bitaille et ah (2007) reported wide ranges of stellar parameters by SED htting: for example. 


protostellar temperature between 3060 and 5013 K and mass between 0.11 and 2.03 Mq. 
Andrews & Williams (2007) observed FT Tau with an angular resolution of 1.7" x 1.1" at 


A = 880 /im using SMA. Although they did not resolve the disk structure, they constrained 
disk properties by model fitting with the data and SEDs. The study of |Guilloteau et ^ 


(2011) using PdBI included this target taken at about 0.5" x 0.3" angular resolution. 


Haro 6-13 has been detected with a continuum flux of 124 ± 13 mJy by Beckwith et ah 
(1990). Schaefer et al.[ (2009) observed the target in CO J = 2 —)■ 1 and J = 1 —)■ 0 as well 
as A = 1.3 mm and 2.7 mm continuum using PdBI. They estimated the protostellar mass 
as 1.0 ± 0.15 Mq, based on the rotational velocity held. They discussed that the difference 
from the mass estimated by the evolutionary track on the HR diagram (~ 0.6 Mq) may be 
due to limited number of channel maps for the CO data model htting. The PdBI data of 
Schaefer et ah (2009) were also used in Guilloteau et al. (2011[). 


HL Tau has been studied very broadly over the last two decades, as it is a bright T 
Tauri star. It has been studied by submillimeter and millimeter imaging of interferometers 
(e.g., Kwon et al. 2011, and references therein), as well as molecular line observations. In 


particular, by a modeling of the viscous accretion disk model to millimeter CARMA data and 
SED data, we found that HL Tau experiences stratihed grain settlement and has an outer 


disk region gravitationally unstable (Kwon et al. 2011). In fact, this paper uses the same data 


set of HL Tau but we employ two disk models and a better temperature distribution. HL Tau 


has also been studied by PdBI (Cuilloteau et al. 2011) and has recently been observed by 


ALMA at three wavelengths (ALMA Band 3, 6, and 7) with an excellent angular resolution 
down to ~ 0.02". 
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3. Observations and Data 

We carried out a survey of protoplanetary disks using CARMA, which is a heterogeneous 
array consisting of six 10.1-m, nine 6.4-m, and eight 3.5-m antennas. For this project the 3.5- 
m antennas were not used for imaging. All targets were observed in A and/or B conhgurations 
(see Table 1^ at A = 1.3 mm and/or 3 mm, which provide high angular resolution better than 
0.3" at A = 1.3 mm. Table provides details on the hnal combined data, their synthesized 
beam sizes and noise levels. 

In addition to the most extended A and B conhgurations, we have observed the proto¬ 
planetary disks in C and D conhgurations, in order to detect large scale structures as well as 
small scales. The highest angular resolution (the smallest synthesized beam size) is ~0.13" 
at A = 1.3 mm toward HL Tan. The angular resolutions for the other targets are a bit less 
(~0.3"). The uv coverage of our data are summarized in Table As shown in the table, 
the minimum uv distance is at most 17 kilo-wavelengths, which means that our large scale 
sensitivity is at least 10 arcseconds. As the target disks are less than a few arcseconds in 
size, there is no resolving-out hux issue in the study, from which interferometric observations 
typically suher for extended targets. 

The largely extended arrays (A and B conhgurations) require stable atmospheric con¬ 
ditions due to long baselines. It is challenging to obtain such good weather conditions even 
at the CARMA site. To verify our data quality, therefore, we have observed an additional 
calibrator (a test calibrator) with a gain calibrator. Since a test calibrator is chosen to be a 
point source nearby our targets, data quality after applying gain solutions can be examined 
by the test calibrator’s shape and hux. All of our A and B conhguration data and some of 
the compact array data have such a test calibrator. Most of our data have built a point or 
point-like image of a test calibrator, which means successful calibration. Here “point-like” 
means that the deconvoluted size of a test calibrator is less than half of a synthesized beam 
in both major and minor axes. A few of the worst data sets have a test calibrator size 
comparable to the synthesized beam and those sets are indicated in Table In addition, 
we examined the hatness of the uv amplitude slopes of the gain calibrators. Since uv data 
are the Fourier transform of an image, a gain calibrator of a point source should have a 
hat amplitude versus uv distance. If it is not achieved, i.e., amplitudes at long uv distance 
regions (small scales, central regions) drop, the target data are likely to be biased, resulting 
in an apparent shallower gradient of density and/or temperature. All of our data sets have 
a hat amplitude slope of gain calibrators after calibration, except a data set taken in B 
conhguration at A = 1.3 mm toward DO Tau, which has about 25% variation in the middle. 
However, it does not seem to ahect the results signihcantly, as the variation is acceptable in 
magnitude and is positioned in the middle of uv coverage. Therefore, we include that data 
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as well with caution, in order to have long baseline data. 


CARMA itself has a special technique for calibration of atmospheric perturbation, called 
CARMA Paired Antenna Calibration System (CPACS) ( Lamb et ah||2009 ). In large antenna 
spacings such as A and B conhgurations, the 3.5 m antennas are paired with boundary 10.4 
m or 6.1 m antennas, which are elements of long baselines. The 3.5 m antennas are always 
observing a calibrator at 30 GHz, while the other antennas are observing science targets and 
calibrators. During the calibration of data, the phase of the 3.5 m antenna gain solution is 
used to correct the phase delay (phase difference) caused by short atmospheric perturbations 
(~4-20 seconds) on the paired antennas. Since visibility phase depends on time delay and 
wavelengths/frequencies, the phase of the gain solution is scaled by the ratio of observation 
frequencies to 30 GHz. The A conhguration data of HL Tan at A = 1.3 mm and the A 
conhguration data of FT Tan at A = 2.7 mm are calibrated using GPAGS. The calibrators 
used for GPAGS are separated from HL Tan and FT Tan by 9° and 13°, respectively, and 
the improvement of calibrated data is about 10-20% in terms of image noise levels and the 
size, flux, and peak intensity of test calibrators. 


We obtained two wavelength data (A = 1.3 mm and 2.7 mm) in order to better constrain 
the dust properties. The properties constrained by multi-frequency data mainly depend on 
absolute flux calibration. To minimize the bias induced by flux calibration uncertainty, a 
good flux calibrator (e.g., Uranus) was used (Table [^. In addition, gain calibrator fluxes of 
all tracks have been compared with each other in the time basis. Rapid variation of gain 
calibrator fluxes (e.g., 50% increase or decrease within a few days) is unrealistic, although 
quasar gain calibrators are intrinsically variable. For comparison, the GARMA flux catalog 
and the SMA calibrator list have been taken into account, as well as the gain calibrator 
values bootstrapped from flux calibrators. In addition, different array-conhguration data 
of the same targets have been compared at common uv distances. Sub/millimeter dust 
emission from T Tauri disks is not variable over a period of a few years, so amplitudes at 
common uv positions should be comparable even in different conhguration arrays. Actually, 
this comparison is a crucial step when combining various array-conhguration data. As a 
result, we presume absolute hux calibration uncertainty of our T Tauri disk data is less than 
the uncertainty assumed in other studies, which use only one or two tracks in an array- 
conhguration for each frequency data: 10% at A = 1.3 mm and 8% at A = 2.7 mm (e.g.. 


Kwon et ah 2011). 


We also considered proper motions of our targets when combining our data, as they 
have been taken over a few year period. However, since proper motions of all targets have 
not yet been well estimated, we decided not to correct them blindly. Instead, we estimated 
the proper motion based on the data themselves. In the case of Haro 6-13, whose peaks show 






an offset in the data set, we matched the peak positions. As a result, we shifted the 3 mm 
data sets by (R.A., Dec.) = (0.0105", 0.055"). In addition, for HL Tan, whose data have the 
highest angular resolution and whose proper motion is well measured, we compensated for 


the proper motion (Kwon et ah 2011). 


MIRIAD (Sault et ah 1995) has been employed to calibrate and map data. Individ¬ 


ual tracks have been calibrated and verified separately then combined in the invert-Fourier 
transform step to make maps. The sensitivity and emphasized scales in the maps depend on 
the visibility data weighting scheme. Natural weighting gives the highest signal-to-noise and 
the worst angular resolution, since visibility data points are sparse at large uv distances and 
visibilities at larger uv distances are noisier due to atmospheric turbulence. On the other 
hand, uniform weighting emphasizes small scales best although the signal-to-noise is worst. 


For intermediate weighting, Briggs introduced a robust parameter (Briggs 1995): 2 of the 


parameter for weighting close to the natural weighting and -2 close to the uniform weighting. 
In order to have small structures reasonably emphasized, the weighting scheme with a robust 
value of 0 has been applied. 


4. MODELING 


We employ two disk models to estimate disk physical parameters: a viscous accretion 


disk model (e.g., Pringle 1981) and a power-law disk model (e.g., Beckwith et al. 1990). 
The outer boundary of the two disk models are described differently. While a sharp outer 
radius is associated with the power-law disk model, a characteristic radius is introduced in 
the viscous accretion disk model, which has a power-law density structure tapered by an 
exponential function. We also consider vertical structures for investigating disk thickness 
and flare degree. For our models, we use cylindrical coordinates, and the disk models are 
axisymmetric. Therefore, the physical properties are expressed in R and coordinates: e.g., 
density structures p{R,z). On the other hand, r indicates a direct distance from the disk 
center, i.e., r'^ = R"^ + z^. The model emission is calculated with a logarithmic radial grid to 
ensure accurate values in the inner region of the disk. 
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4.1. Power-law Disk Model 


The power-law disk model has a power-law density distribution with a sharp edg(J3 


p{R,z) = p{R,^) ey.Y>[-{z/H{R)f] 

= exp[-(z/i7(/2))2], (1) 

where H(R) is a scale height along radius set by a power-law function: 

= ( 2 ) 

Here Hq = {2kR^TQ/GM^fh)^'^ is the hydrostatic equilibrium case between the vertical 
component of the protostellar gravity and the local pressure of T(i?0)0) = Tq. We ignore 
self gravity of the disk mass and assume that dust temperature and gas kinetic temperature 
are the same. In the Hq expression, k is the Boltzmann’s constant, G is the gravitational 
constant, M* is the central protostellar mass, and rfi is the mean molecular mass: fh = 
pmjj = 2AmH- Note that disks in the hydrostatic equilibrium and in the optically thin 
condition have 6* = 1 and h = 3/2 — q/2 zz 1.25, where T(i?, 0) oc 

As the surface density is obtained by the integration of the density distribution in z, 
the surface density is also a power law with a power index oi p = s — h\ 

/ OO 

p(i?, z)dz 

■OO 

= p(R,0)V^H{H) 

= p„^hH„(Ay, ( 3 ) 

In this model we also have a sharp inner radius (i?m) in addition to the outer radius (Rout)- 
Therefore, po is expressed with a disk total mass {Mdisk), 


^disk 


Po 


f* Rout /*oo 


p{R, z)dz 2 'kR dR 


' Ri . 


Po 


M, 


271^-^RlbtHo 
2 — p 

2 — p 


R 


out 


Rn 


2-p 


Rin\^-P' 

Ro 


disk 


271^-^RlbtHo R 


R. 


out 


2-p 


Rin\^-P 


Rn 


-1 


( 4 ) 


'^We use s for a volume density distribution and p for a surface density in order to minimize confusion 
with other studies. Note that this notation is different from 


Kwon et al. 


( 20111 . 
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Although we consider a power-law temperature distribution for the disk scale height, we 
use a temperature distribution calculated in a given disk structure by the RADMC 3D code 
( Dullemond fc Dominik||2004 ) for more realistic modeling with a temperature gradient along 
the vertical direction. However, we did not update the scale height based on the calculated 
temperature distribution. This decoupling of the disk structure from the temperature dis¬ 
tribution is justified, because dust grains settle to the mid-plane rather than mix with gas 
as they grow (e.g., Andrews et al. 2009). 


The RADMC 3D code provides two dimensional temperature distributions, but it re¬ 
quests a large number of photons for a smoothly varying distribution particularly in the case 
of thin disks with a high optical depth. Therefore, instead of simply using the tempera¬ 
ture distributions out of the RADMC 3D code, we make a functional form weighting two 
functions that fit to the mid-plane {Tm) and surface temperature distributions (T<j), 


r(R, z) = WTm{R, 0) + (1 - W)Ts{r). (5) 

Based on qualitative comparisons, we chose the weighting function W = exp[—( 2 :/if(R))^]. 
The dust opacity spectral index f3 is one of the most interesting parameters in this 


study, which represents dust grain sizes (Draine 2006): Ky = kq{i'/i'qY. Indeed, in order to 
constrain /3, we fit the two wavelength data simultaneously. Therefore, we assume kq = 0.01 
cm^ g“i at 230 GHz and (3 is determined between A = 1.3 mm and 2.7 mm. The kq is based 


on the case of ice mantle grains following a MRN size distribution (Mathis et al. 1977) 


m 


Ossenkopf & Henning (1994). In addition, we assume a gas-to-dust mass ratio of 100. 


The assumed kq corresponds to the case of Ky = 0.1(i//1200 GHz)^ with /3 = 1.39 (e.g., 
Hildebrand||1983 Beckwith &: Sargent||1991 ). The uncertainty of kq is large, around factor of 
two ( Ossenkopf fc Henning||1994 ). For the opacity as an input of RADMC when determining 
disk temperature distributions we simply consider the combination of an opacity curve at 
infrared wavelengths in the case of silicate-carbon grains in the MRN size distribution and 
a power-law at millimeter wavelengths. 


4.2. Viscous Accretion Disk Model 


In addition to the power-law disk model, we utilize the viscous accretion disk model for 
fitting our data. The viscous accretion disk model has a density distribution of a power-law 


tapered by an exponential function (e.g., Pringle 1981), 


Rr 



r/ R \ 2 -S +/11 


r f ^ Y] 

1 exp 


exp 

[ \H{R)) \ 


( 6 ) 
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The viscous accretion disk model is expressed as the surface density distribution (e.g., An¬ 


drews et ah 2009), 


E(fl) = 


p(i?, z)dz 


= poy/nbtHo ^ — 


Rc\^ / R \ 


Rn / V Rr / 


exp 


/ R \ 2-S+/X 


= exp 


RcJ 


R \2-7] 


RcJ 


( 7 ) 


where Rc is a characteristic radius. Note that Sq = Po\/ttH{Rc) and 'y = s — h. Like 
the power-law disk model, we decouple the scale height distribution from the temperature 
distribution, and the density distribution is expressed with the total disk mass. 


J^disk 


Po 


POO roo 


p{R, z)dz 2 ttR dR 


Rin 

2n^-^RlHiRc) 

n c \ 


^ Rin \ 

2-7- 

2-7 


\RcJ 

- 

2-7 

T\/T.. . ' 

- 

f Rin ^ 

2-.. 


V Rc t 

J 


( 8 ) 


Everything except the density distribution is set in the same way to the power-law disk 
model. Note that the characteristic radius (Rc) is introduced instead of the outer radius 
( Rout ) • 


4.3. Model Parameters 

There are 7 free parameters for each disk model: s, /3, Mdisk, Rim Rout (or i?c), di, and 
PA. In addition, there are two parameters related to the disk scale height, h and bt, which 
are weakly constrained by the observations and will be handled as “controlled parameters” 
in the htting. The 6^ is an inclination angle of a disk (a face-on disk with 6i = 0°) and PA 
is a position angle of a disk measured eastward from the north. The Rout and Rc are for 
the cases of power-law disk and viscous accretion disk models, respectively. For the viscous 
accretion disk model, an outer disk radius is not employed. However, the modeling disks 
are cut off at 5 Rc due to limitation in image size and for simplicity in integration. Some 
modeling studies of circumstellar disks have presumed sublimation radii of dust grains (e.g., 
, but it is possible that inner radii {Pin) are determined by other effects 
such as an unresolved binary companion, formed planets, etc. Therefore, we leave the inner 
radii as a free parameter, but htted values can be model-dependent and are likely unreliable 
due to systematic uncertainties in the observation, primarily the resolution. 


Andrews et al.||2009) 
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The sign of inclinations are determined by our spectral line data in CO, which were 
included in our observations. We only detected CO in data from the C and D conhgurations, 
which traces the bipolar outflows as well as a disk structure. This allowed us to estimate 
the disk inclination, at least the sign. Cl Tan and FT Tan, which we have not detected any 
CO, are set to have a plus sign of inclination. All the inclination signs are consistent with 
previous studies (e.g., Guilloteau et al. 2011). 


The two controlled parameters, bt and h, are used to investigate the disk vertical struc¬ 
ture. As shown in equation bt and h are disk scale heights with respect to the hydrostatic 
equilibrium cases and disk flare indexes. In order to hgure out whether preferred disks 
are thin or thick and more or less flared, we carry out 25 separate models for each target: 
bt = {0.6, 0.8,1.0,1.2,1.4} and h = {0.95,1.10,1.25,1.40,1.55}. We adopted the 25 cases 
instead of having them as two additional free parameters for maximizing the efficiency of 
our modeling without losing our study goals. 


4.4. Model Fitting Procedures 


First, a temperature distribution is obtained by the RADMC 3D code in a given disk 
structure and fitted as described in Section 4.1 Then, a disk model is built in two dimensional 
logarithmic grids by solving the radiative transfer equation numerically along the line of 
sight. Note that the disk inclination is applied when constructing the disk model, so the disk 
model must be two dimensional and the radiative transfer equation is numerically integrated 
along the third axis, line of sight. No optically thin nor Rayleigh-Jeans approximations 
are assumed. Afterward, disk images are made by interpolations of the disk models in two 
dimensional linear pixels. The center positions and position angles (PAs) of the disks are 
applied when making the disk images. 


The disk images are multiplied by three different normalized primary beams, which 
correspond to three different types of baselines in CARMA: 10.4 m-10.4 m antennas, 10.4 
m-6.1 m antennas, and 6.1 m-6.1 m antennas. These three primary beam corrected images 
are Fourier-transformed into model visibility maps. Then, the actual model visibilities are 
sampled in the uv coverage of observational data by bi-linear interpolation of the model 
visibility maps. Image pixel sizes and image sizes are selected to obtain reasonable pixel 
sizes and uv coverage of visibility maps. For example, an image pixel size of 9 AU and an 
image size of 1024 pixel by 1024 pixel are used, which provide uv coverage up to 1600 k\ {uv 
range: -1604 to -1-1604 kX) and a visibility map pixel of about 3.1 kX, based on the Nyquist 
theorem: for example, the uv range A^v = ^/26image, where Stmage is the image pixel size. 
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Model fitting is done by comparing observational data with the model data sampled 


along the observational uv coverage in Bayesian inference (e.g., Gilks et ah 1996; MacKay 


2003). Bayesian inference allows us to obtain the probability distribution of disk properties 


(m) with given data D and a given disk model H: P{m \ D, H), 


P{m \D,H) = 


P{D I m,H)P{m \ H) 
P{D I H) ' 


(9) 


Here P{m \ D,H) is called the posterior probability, P{D \ m,H) the likelihood, P{m \ H) 
the prior, and P{D \ H) the evidence. The evidence is also called the marginal likelihood, 
which is used to compare multiple models, since P{D \ H) = f dmP{D \ m, H)P{m \ 
H) ^ E P{D I m, H)P{m \ H). Gaussian functions are employed for the likelihood, 
since noise of interferometric data is a normal distribution, ll.exp[-{Di- The¬ 
oretical uncertainties of interferometric visibility points are typically estimated as cr = 


{2k/Arf)^^TiT2/2tint^i^ (e-g-, Thompson et ah 2001), where k is Boltzmann’s constant, A is 
antenna area, rj is the efficiency of the antenna surface and correlator quantization, Ti and 
T 2 are system temperatures of antennas corresponding to the baseline, Au is bandwidth, 
and tint is integration time. However, the uncertainties do not represent decorrelation of 
atmospheric turbulence, which depends on uv distance. Therefore, instead we use the stan¬ 
dard deviation of the imaginary components of self-calibrated gain calibrator data, which 


are generally a few times larger than the theoretical estimate (Kwon et ah 2011). Individual 


wide bandwidth windows {Au ~ 500 MHz) are averaged and considered a visibility point. 
The central channel frequency of the window is used for computing its uv distance in units of 
wavelength. As the data are real numbers in image space, i.e., the uv coverage is symmetric 
at the phase center, the symmetric data points are added to the sample. In addition, the 
real and imaginary components of complex visibility data are regarded as independent data. 

Searching distributions of model parameters by the Metropolis Hastings method has two 
parts. The hrst part is finding a convergence region of parameters. For a better performance, 
we start with 8 initial conditions, which are randomly generated. After running for a while, 
the initial conditions converge into a parameter set (normally two sets with different signs 
of inclination angles). Annealing is also used to accelerate the convergence. We start a fresh 


run from the convergent parameter set of the expected inclination sign (refer to Section 4.3) 
with a smaller standard deviation of the proposal Gaussian functions, for each combination 
of bt and h, 25 in total for each target. The smaller standard deviation of the proposal 
functions are set to have an acceptance rate of about 10-50 %. The GASDEV subroutine 
of the Numerical Recipes using the Box-Muller method is used to produce random numbers 


in the Gaussian distribution (Press et al. 1996). 
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5. Results and Discussion 


5.1. Combined Data and Images 


Figure [T] shows A = 1.3 mm continuum images of the six protoplanetary disks; a physical 
scale bar of 100 pc and the data angular resolution (synthesized beam) are marked. Note 
that the three disks in the bottom panels (Cl Tan, DL Tan, and HL Tan) are extended, 
while the three disks in the upper panels (Haro 6-13, DO Tan, and FT Tan) are relatively 
compact. In addition. Cl Tan has a weak non-axisymmetric feature north-north-east to 
south-south-west. 


The continuum noise levels that we achieved after combining all configuration data 
are listed in Table with disk total fluxes. In addition. Table has disk center positions 
estimated by a Gaussian £t (the IMFIT task of MIRIAD). Table has phase center coor¬ 
dinates, so one can simply compute the actual disk centers: e.g., HL Tan, R.A. (J2000) = 
04:31:38.418, Dec. (J2000) = +18:13:57.37. 

Figure 1^ presents observational visibilities of six targets averaged in annulus with the 
best fitting model overlaid: the solid and dashed lines indicating accretion disk and power- 
law disk models, respectively. The dust opacity spectral indexes at the bottom of each panel 
have been calculated in the optically thin assumption and the Rayleigh-Jeans approximation: 
/3 ^ \og{Firam/F 3 mm)/^og{uimm/j^ 3 mm) “ 2. The /3 values from our modeling without such 
assumption and approximation are addressed later. Open squares and open triangles mark 
A = 1.3 mm and 2.7 mm data, respectively. The error bars are the standard error of the 
mean: where af is the variance of a visibility data point in a annulus calculated 

from its Tgyg and N is the total number of data in the annulus. Absolute flux calibration 
uncertainty is not applied for the plots. As shown, the two models are htting observational 
data well. Note that since the viscous accretion disk model does not have an outer radius, 
there are no bump features at around 200 — 300 k\, corresponding to about 1". 


Although the f3 values of the figure are simply calculated and not from our modeling, 
it is worthy of noting the weak relative variations along uv distance, which can be used for 
studying grain size distributions. In general, recent studies reported that protoplanetary 
disks have larger grains in the inner regions (e.g., Guilloteau et al. 2011 Perez et al.|2012 


Trotta et al. 2013). Observations with higher angular resolution and sensitivity at multiple 


wavelengths are required to investigate the spatial distributions of dust grains. We do not 
intend to study the grain size distributions in this paper and simply assume a constant (3. 












15 


5.2. Fitting Results 


We have carried out two flared disk models: viscous accretion disk and power-law disk 
models. As described in Section we have 7 free parameters and two controlled parameters 
for the vertical structures (5 disk thickness bt and 5 flare index h values). The fitting results 
are presented in Tableand Figures]^ andThe 25 cases of the controlled parameters are 
combined in order to see the overall parameter ranges, as they have been sampled with the 
same Gaussian proposal function (Section |4.4[ ). 

We did not use any priors obtained from the images of the targets except the center 
positions and the inclination signs. We tested the center positions determined by a Gaussian 
fit of the IMFIT task in MIRIAD and the ones manually determined, based on the intensity 
peak positions. The offset values by the simple Gaussian fit have been chosen (Table |^, as 
they give higher posteriors. 

The fitting results are listed in Table and the posterior distributions of parameters 
in individual targets are presented in Figure and The table values are the means and 
standard deviations of individual parameters in the marginal probability distributions. The 
solid lines of Figure to indicate normal distributions. As shown, most parameters have 
a normal distribution posterior, and hence well defined values and deviation. 


The exception is the inner radius. The inner radius was fit over a range from 1.3 to 9 AU 
because none of the disks showed an inner hole or central depression of emission at our best 
18 AU resolution. The resulting fits show weak dependence on inner radius. The statistical 
uncertainties in Table do not acconnt for the systematics of increased atmospheric de¬ 
coherence on the longest baseline which are biased toward blurring the central disk emission. 
The inner radius for HL Tau, the source with the best resolution and signal-to-noise, may 
represent a flattening of the surface density distribution within 10 AU rather than a central 
hole (i.e. ALMA SV data). 


The volume density distribution index (s) of the power-law disk model is between ~2.0 
to ~2.8, which corresponds to a surface density index (p) between ~0.4 and ~1.8 when the 
best flare index values are adopted {p = s — hm. Section]^. On the other hand, the s values 
of the accretion disk model are in a range of 0.7 to 2.2, which are converted into 7 expressing 
the surface density distribution (7 = s — h) of -0.2 to 1.1, as shown in Table Such a large 


range of 7 has been reported by previous studies (Isella et ah 2009 Andrews et ah 2009). 


The dust opacity spectral index (3 appears to be constrained very closely in the two 
models. This is expected because (3 is sensitive to the difference of flux densities at two 
wavelengths. (3 of our sample disks ranges from -0.15 to 0.70 and half of them are around 


0, which can be interpreted as grain growth. 13 of younger Glass 0 YSOs is around 1 (Kwon 












et al. 2009|), so grains may grow with object evolution. This is still valid for DL Tan, FT 

Note that the possible error of [3 due 

at A = 2.7 mm; 


Tan and HL Tan, since they have (3 smaller than 1. 
to the absolute flux calibration uncertainties (10% at A = 1.3 mm and 
refer to Section]^ are about 0.25. 

Disk masses are 0.01-0.10 Mq. Compact disks (DO Tan, FT Tan, Haro 6-13) tend to 
have smaller masses. However, disk masses are distributed in a large range, among which 
HL Tan has the largest mass. Note that the disk masses have a large uncertainty mainly 
because kq has a factor of two uncertainty (e.g., Ossenkopf & Henning 1994). In addition. 


the stellar luminosity affects the disk mass, as it sets the dust temperature. However, disk 
masses are not so sensitive to the stellar luminosity due to the low power index and the high 
optical depth along radius. Note that the derived temperature distributions for our modeling 
are calculated by the Monte-Carlo radiative transfer code of RADMC 3D. 

The disk sizes (characteristic and outer radii) are well constrained in both models: 
Rc ~ 30-120 AU and Rout ~ 60-200 AU. The characteristic radius of the accretion disk 
model is where the disk density distribution changes from the power-law dominant region 
to the exponential dominant region (Equation]^. Comparing the surface density equation 
of our accretion disk model with a similarity solution of viscous accretion disks and con¬ 
sidering the mass flow equation (e.g.. 


Andrews et al. 

2009 

Hartmann et al. 

1998 

jPringle 


1981), the characteristic radius can express the transitional radius 7?^, where the bulk flow 
direction changes from inward to outward to conserve angular momentum (e.g.. Equation 


A9 in Andrews et al. 2009): 


Rf = Rr 


1/(2-7) 


■L2(2-7). 

The Rt of our sample disks are in a range of about 15 to 60 AU, as listed in Table 


( 10 ) 


The two geometrical parameters, inclination (0j) and position angle (PA), are the two 
best constrained. The values are consistent with previous high angular resolution observa¬ 
tions. In addition, they are consistent in both models, so model-independent. We originally 
did not limit the inclination sign based on the known bipolar outflow directions of our targets 
in the modeling. However, for the final modeling results, we have chosen the inclination sign 


in agreement with the bipolar outflow direction and/or previous studies (e.g., Guilloteau 


et al. 2011) 


Five of our disks were also observed and fit to power-law and accretion disk models by 


Guilloteau et al. (2011). Disk inclination and position angles agree fairly well (accounting for 
the difference in definition of position angle) but other disk parameters show only broadly 
similar trends, without detailed agreement. For HL Tan and Haro 6-13, our angular resolu¬ 
tion, which is a factor of 3-4 higher, is the likely reason for the difference. For Cl Tan, where 
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the resolutions of the observations are nearly the same, the htted disk parameters, Re, 7 , 
Rout, and p, are similar, but generally outside of the statistical uncertainties by more than 
2a. Guilloteau et ah (2011) assume a simple R~^'^ power-law for the radial temperature 
variation which may be responsible in part for the difference. However, it is likely that both 
analysis are more limited by the sensitivity and angular resolution of the observations than 
is reflected in the statistical uncertainties. This is a cautionary note for all disk modeling to 
date. 


5.3. Model Comparison 


First, we compared posteriors of the two models. Bayesian inference provides a way 
for model comparison; models can be compared by evidence, which is the integration of 
likelihood times prior over the whole possible parameter space of a model (e.g., MacKay 


2003). We set the prior as well as the number of free parameters the same for the two 
models (uniform priors over the same ranges), so we just need to integrate likelihood. The 
likelihood integration can be estimated by the best likelihood times the posterior accessible 
volume: 


K = 


P{D I Ha) 
P{D I Hp) 


J dmP{D 

m, HA)P{'m 

\Ha) 

f dmP{D 

m, Hp)P{m 

\Hp) 


P(D 

r«s_; 


HA)(^{m\D)A 

P{D\ 

P(D 

r«s_; 

^^best-! HP 

‘^best'i 

Hp)a(^m\D)p 

P{D\ 

latest, Hp)(^(m\D)p 


( 11 ) 


where the subscript A and P indicate the accretion disk model and the power-law disk 
model, respectively. The posterior accessible volume is estimated by cr(m|D) = IfcTj, where exj 
is posterior widths of s, j3, M^isk, Pin, Rout (or Pc), 6i, and PA. 

The K values of our targets are shown in Table with the evidence of the two models. 
As shown in the table, all the targets, except DL Tan, have a positive value in In(iF). 
Particularly, HL Tan, whose data have the highest angular resolution in our sample, has 
a signihcantly larger positive value. A scale for interpretation of these values has been 
proposed: K > 3 indicates a considerable preference for the model of the numerator (here 
the accretion disk model) and K < 1/3 for the denominator model (here the power-law 
disk model) (e.g., Kass & Raftery 1995[ ). Therefore, we conclude that overall the accretion 
disk model is preferred. It is noteworthy that the preference of the accretion disk model is 
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quantitative and purely based on continuum data. Hughes et al.| (2008) qualitatively argued 
that the accretion disk model with a tapered outer region instead of the power-law disk model 
with a sharp outer radius is preferred, as it explains more extended disk features detected in 
gas tracers. A previous work on two protoplanetary disks reported slightly smaller reduced 
values for the accretion disk model than the power-law disk model (Isella et al. 2010). 


In contrast to the overall trend, DL Tau prefers the power-law disk model, which may 
be implying a difference in evolutionary state or physical conditions in the disk. In addition, 
the model preference is not signihcant toward Cl Tau and Haro 6-13. It may be due to the 
non-axisymmetric feature in Cl Tau and the small disk size of Haro 6-13. Note that Haro 
6-13 has the smallest disk in our sample, slightly larger than our angular resolution. 


Guilloteau et al. (2011[) also reported no overall preference between the two disk models 


in observed disks with sizes comparable to the angular resolution. However, they found 
that HL Tau prefers a power-law disk model, while DL Tau prefers an accretion disk model, 
unlike our results. Although their data have slightly lower angular resolution than ours, they 
took the radial dependence of {3 into account. We have only a constant (3, but we employ 
a better temperature distribution and a vertical structure. The discrepancy simply presents 
the difficulties of model comparison using the current data sets. 


In fact, disk model comparison studies may require higher sensitivity and resolution 
data for a larger sample, as well as detailed modeling, in order to achieve a reliable result. 
We do not attempt further discussion of the model comparison results in this paper. Instead, 
we further investigate the disk properties obtained only from the accretion disk model in the 
following. 


5.4. Residual Maps 


Figure [^presents the residual maps for the best accretion disk models with the A = 1.3 mm 
and 2.7 mm continuum images: observational images, best models, and residual maps from 
the left. The residual maps do not have any signihcant features except FT Tau and HL Tau. 
The scattered peaks of FT Tau at A = 2.7 mm are mainly due to the noisy data of the longest 
baselines taken in A-con£guration. The residual features of HL Tau have been discussed in 


Kwon et al. (2011): the A = 2.7 mm features may be due to free-free emission contamination 


and the A = 1.3 mm features are hints of substructures now seen in the ALMA SV data. 
The low levels of the residual indicate excellence of our model htting. 


Although the residual maps illustrate that the data are consistent with our models, the 
well-resolved map of Cl Tau has a weak substructure that is seen at both wavelengths; a non- 
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axisymmetric feature is present from north-north-east to south-south-west at A = 1.3 mm 
and from north to south at A = 2.7 mm. Although the feature is weak, the similar pattern 
at both wavelength data sets provides a better hdelity for a non-axisymmetric substructure. 
DL Tau and HL Tau do not show a distinct feature for a substructure, but the extended 
disks include a possible region for the gravitational instability. We discuss the possibility 
of substructure development on these disks with the Toomre Q parameter later. The HL 
Tau data and images have been presented in Kwon et ah (2011), but we have employed a 
slightly different accretion disk model with an independent disk flareness parameter and a 
more realistic temperature distribution in this paper and carried out modeling identically 
for all the six targets. 


5.5. Correlations among Properties 


We calculated correlations between the constrained physical properties using 


R’xy 


T.{xi - x){yi - y) 


^jT.{xi - xY^jT.{yi - yf 


( 12 ) 


where x and y are means of {x,} and {^i}, respectively. The correlation coefficients are 
presented in Table Relatively strong correlations (> 0.5) are emphasized with a bold font 
and plotted in Figure]^ Although these correlations are obtained in our small sample of six 
targets, a few relationships may provide meaningful hints on disk properties and evolution. 
In addition, note that the four correlations are the strongest (> 0.5) even without HL Tau 
data points that look dominant in the relationships. 

As we obtained multiple wavelength continuum data and £t the data simultaneously 
with a single model, we can estimate the dust opacity spectral index (/?). Overall, (3 of 
our sample protoplanetary disks are 0.6 or less, which indicates that grains have grown 


larger than the younger protostellar systems (Kwon et ah 2009), as mentioned in Section 
|5.2[ Previous studies using high angular resolution data over a broad wavelength region 
(submillimeter to centimeter) have reported radial dependence of (3 (e.g., Perez et ah 2012). 


Our model htting does not show evidence for such a radial dependence in (3. In specihc, 
our sources with the most complete data out to large uv distance at 1.3 and A = 2.7 mm, 
DL Tau, Haro 6-13, and HL Tau, have no radially symmetric residuals in Figure 5. This is 


in conflict with Guilloteau et ah (2011) which finds signihcant radial variation in [3 which 


should have been seen in our data. The (3 versus radius for our best resolved disk, HL Tau, 
is not given in Guilloteau et ah (2011). Higher resolution and high sensitivity observations 


over a larger range of wavelengths are needed to solve this discrepancy. 
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We found a correlation between f3 and the disk masses: disks with a smaller /3 are less 
massive. Small indicates large grains. Therefore, this relationship implies that less massive 
disks have a larger fraction of large grains. One caveat on this correlation is that {3 close to 
0 may result in an underestimate of the disk mass because k could be significantly different 
from the assumed value. As an extreme, /d = 0 is consistent with any size dust grains larger 
than roughly an observational wavelength - pebbles, rocks, etc. This relationship has also 


been found in circumstellar disk around the more massive young Herbig AeBe stars (Sandell 


et al. 2011 ). 


We also found that the mid-plane density (s) is anti-correlated with the opacity spectral 
index (/d) and the disk mass In other words, protoplanetary disks with smaller /3 

and of a smaller mass have a steeper mid-plane density gradient. This indicates that disks 
with large grains tend to be detected with lower mass and to have relatively larger mass in 
the inner regions. The increase of the mid-plane density gradient suggests that large grains 
move radially as well as settle down in the mid-plane. However, (3 does not have a strong 
relationship with the flare index h, which implies that grain growth can occur independently 
of disk flatness and grain settlement. It is noteworthy that there are other mechanisms of 
grain growth and disk flatness such as grain drift and disk evaporation. On the other hand, 
grain settlement is largely affected by turbulence, which can be studied by disk gas tracers. 
Interestingly, the flare index does have a relationship with disk size as shown in Figure and 
described in Section 15.61 


However, we did not obtain anti/correlations supporting the ideas that the equivalence 
width of small grains (EW) and the spectral index between 13 and 31 /rm ( 77 . 13 - 31 ) decrease 
when grains grow and disks flatten. Conversely, we obtain the opposite trends (see Table 
[^. For example, (3 is anti-correlated with EW, which implies that disks with large grains (a 
smaller (3) have more fine grains (a larger EW). In the sense that our millimeter wavelength 
data are sensitive to large grains and infrared data are sensitive to hne grains, it can be 
interpreted with a more dusty disk. In addition, note that disks can be detected as a thick 
disk at short infrared wavelengths and as a thin disk at long millimeter wavelengths: stratihed 
grain settlement (Kwon et al. 2011). 


5.6. Disk Vertical Structures 

Figure shows the posterior distributions of the accretion disk models in the disk 
thickness (bt) and the flareness index (h), normalized by individual peak values. As shown, 
the three extended disks (Cl Tan, DL Tau, and HL Tau) have higher posteriors in thinner 
and less flared models. In contrast, DO Tau appears to be better £t by thick and more flared 
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models. In the cases of FT Tan and Haro 6-13, distinct preferences are not present. We 
argue that the sensitivity and angular resolution of the data are not good enough for these 
objects. Note that they are the two faintest targets and compact disks of our sample. 


While our modeling results of millimeter wavelength data prefer thinner and less flared 
disks for Cl Tan, DL Tan, and HL Tan, a thick disk model can be required for shorter mid- 
IR observations. For example, HL Tan needs a thick disk for its mid-IR fluxes, as shown 
Kwon et ah (2011). The discrepancy can be understood by stratihed grain settlement: 


m 


large grains, which millimeter data are sensitive to, have been settled down in the mid-plane, 
while hue grains are not. 


The dust opacity spectral indexes [3 do not have a clear (anti-)correlation with the disk 
thickness/flareness nor the disk sizes, while larger disks are thinner and less flared. Therefore, 
it may imply that large grain settlement is not the only mechanism for thinner and less-flared 
disks. This is consistent with the fact that protoplanetary disks are evaporated and get less 
flared and thinner with time. 


Furlan et ah (2009) argued that protoplanetary disks with grain settlement have small 


?7 ,i3 _3 i and EW(10 /xm) values. Note that EW(10 /im) indicates the amount of small (< 5 /im) 
silicate grains and ni 3 _ 3 i represents the flareness of disks by intervening of protostellar emis¬ 
sion. Unlike their argument, we did not hnd a strong relationship between both properties 
and our flare index (h). Probably, the short wavelength properties of EW(10 /rm) and ni 3 _ 3 i 
are sensitive to the warm inner and/or surface disk region and relatively small grains, while 
our millimeter data are sensitive to the cold outer and/or deep disk region and large grains. 
Also, inclinations of objects may be an effective parameter, as the light path so the optical 
depth along line-of-sight depends on it. In contrast, we obtained a strong correlation between 
s and EW(10 /im) and a strong anti-correlation between (3 and ni 3 _ 3 i. 


5.7. Disk Substructures 


As discussed in Kwon et al. ( 2011| ), HL Tan has a region between 50 and 100 AU, 
which is gravitationally unstable and could fragment. The prediction of substructures in HL 
Tan has been proved by the ALMA SV data with the unprecedented angular resolution up 
to 0.02", which resolves at least 8 gaps (and rings). Particularly, the data present a wide 
gap between 50 and 100 AU. Indeed, our CARMA data show a hint of the wide ring, as 
shown in Figure In addition, the positive and negative residuals nicely line up with the 
rings and gaps. For the comparison of the CARMA residuals with the ALMA image, we 
shifted CARMA data by the proper motion. The HL Tan center positions are (R.A., Dec. 
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in J2000) = (04:31:38.418, +18:13:57.37) in the CARMA observations on 2009 Jan. 31 and 


(04:31:38.42545, +18:13:57.242) in the ALMA observations on 2014 Oct. 30 (Partnership 
et ah 2015), so the proper motion is estimated as vr,a. = 18 mas year“^ and Vdeci. = —22 


mas year“^ during the 5.75 year span. In addition, the temperature distributions of our 
modeling are consistent with the brightness temperature estimated by the ALMA multi¬ 
band SV data. Figure shows examples of our model temperature distributions in color 
solid and dashed lines. In each case, the upper straight line is for the surface temperature 
and the bottom curve or line is for the mid-plane temperature. As described in Section 
the temperatures were obtained by htting RADMC 3D results. The black line is the outline 
of the brightness temperature estimated from the ALMA multi-band SV data (Figure 3.d of 
Partnership et ah 2015). Note that the brightness temperature is between the surface and 


mid-plane temperatures of our modeling, which means that our temperature distributions 
are reasonable. 

Like HL Tan, whose CARMA data show a hint and whose ALMA data provide the 
detailed substructures. Cl Tan and DL Tan images also show a possible non-axisymmetric 
feature. We calculated the Toomre Q parameter to check if the three extended targets have 


gravitationally unstable regions. The Q parameter is dehned as (Toomre 1964), 


ttGE’ 


(13) 


where c* is the isothermal sound speed, D is the orbital angular velocity (D = 

G is the gravitational constant, and S is the surface density, and regions of Q < 1.5 are 
gravitationally unstable. As we constrained disk temperature and density distributions, 
the Q parameter can be calculated assuming Keplerian rotation. Note that we utilize our 
mid-plane temperature distributions for c*, instead of a constant temperature: = kT'^. 


Figure [TO] shows the Q parameter along radius in the accretion disk models. As shown 
in the figure. Cl Tan and DL Tan have the smallest Q minima with HL Tan, but they do 
not have a region with Q < 1.5. However, it is possible that they could have undetectable 


mass in large grains, as discussed for the correlation between j3 and disk mass in Section 5.5 


in particular, /3 of Cl Tan is close to 0. On the other hand, DL Tan has a region with the Q 
values very close to 1.5 around 100 AU in radius, although the sensitivity of our data does 
not provide any distinct non-axisymmetric features. 

Non-axisymmetric substructures are a possible place for protoplanet formation devel¬ 
oped by gravitational instability. Also, they allow us to study the accretion mechanism 
at the late stage of star formation. A trailing spiral structure built by gravitational insta¬ 


bility is an efficient mechanism to transport angular momentum outward (e.g., Hartmann 


2001). Therefore, if we could obtain the kinematics of such a feature along, we can study 
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the accretion mechanism. As the SV data proved, high angnlar resolution and sensitivity 
observations of ALMA will clearly show whether or not the extended disks of our sample 
have a substructure and enable us to study the accretion mechanism. 

While the accretion mechanism of the gravitational instability generally appears in the 
outer disk region, the most successful mechanism to explain angular momentum transport 


in the accretion disks is the magnetorotational instability (MRI) (e.g., Balbus 2003). The 
key aspects MRI requires are magnetic helds coupled with material to act as a tension and 
disk rotational velocities decreasing outward. As with the trailing spiral arms, in which the 
gravity of the mass excess acts as a tension, a faster rotating inner region is dragged by 
the slower rotating outer region connected by the magnetic held tension. Therefore, angular 
momentum is transported outward. The rotation decreasing outward is a general property 
of circumstellar disks, whose velocity appears Keplerian. For the former aspect of magnetic 


helds coupled with disk material, the key property of disks is ionization fraction. Gammie 


(1996) tackled which parts of disks are coupled with magnetic helds and suggested a layered 


accretion disk. While the inner region (< 0.1 AU) is collisionally ionized, the outer region 
is layered with accretionally active (ionized enough by cosmic rays) and dead zones (cold 
mid-plane). 

The viscous accretion disk model provides a means to investigate the physical origin of 
viscosity. Viscosity of disks is normally parameterized as z/ = aCgH, where a is a dimen¬ 
sionless parameter, Cg is a sound speed, and if is a scale height (e.g., Shakura & Sunyaev 


1973). On the other hand, since the mass accretion rate is related to the viscosity in a 


structure-determined disk (e.g.. Equation AlO in Andrews et ah 2009), we can estimate the 
a parameter of our disk models htting data, assuming an accretion rate. For this calculation. 


we use the mass accretion rate estimated by Robitaille et al. (2007) using SED htting (Table 
and the equation derived by Andrews et al. (2009|): 


a 


2Rl r Ry 1 


3(2-7)MrfVRc^ CgH 


(14) 


Although the accretion rates determined by UV and optical spectroscopy are more accurate 


(Robitaille et al. 2007), we adopted the SED htting results in order to have as many target 


values as possible. Note that the values estimated by UV and optical spectroscopy are 
available for Cl Tan, DL Tan and DO Tan and are close to the minima of the SED htting 
results in the parentheses of Table The a is a function of radius but we listed only the 
values at 10 AU and 100 AU. The a values at 10 and 100 AU of our sample disks are 
all in the range that MRI can be the viscosity origin for: 0.005-0.6 (Balbus 2003). We 


assumed the accretion rate of Haro 6-13, which has not been reported, as a typical value of 
1.0 X 10“^ Mq year“^. The investigation of the viscosity origin is beyond the scope of this 
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study, so further attempts to understand the exception are not made. 


6. Conclusion 


We have carried out a T Tauri disk survey using CARMA, which provides excellent 
image hdelity and angular resolution. We have acquired multi-wavelength (A = 1.3 mm and 
2.7 mm) and multi-conhguration (A, B, C, and/or D) data up to angular resolution of 0.13" 
toward 6 targets: Cl Tan, DL Tan, DO Tan, FT Tau, Haro 6-13, and HL Tan. Using visibility 
modeling with Bayesian inference, we obtained disk properties of the flared viscous accretion 
disk model such as density distribution, dust opacity spectral index, disk mass, disk inner 
and characteristic radii, inclination, position angle, and disk vertical structures. In addition, 
we examined anti/correlations between the properties. The power-law disk model was also 
employed for investigating which disk model is preferred. 

1. We found that the accretion disk model is preferred overall, with the exception 
of DL Tau, which prefers the power-law disk model. However, the understanding of the 
discrepancy is not clear at the moment. ALMA observations with a better sensitivity and 
angular resolution toward a large number of protoplanetary disks are required for such model 
comparison and trend investigation further. 

2. We obtained correlations between the properties we constrained using the flared 
accretion disk model. Particularly, we found that disks with a steeper mid-plane density 
gradient have smaller /3 (large grains) and are less massive. This suggests that grains grow 
and radially move and implies that a less massive disk has a larger fraction of large grains. 
In addition, we found that extended disks tend to be less flared. 


3. We detected a non-axisymmetric feature toward Cl Tau with a similar pattern at 
both millimeter wavelengths. However, the Toomre Q parameter criterion does not support 
a gravitationally unstable structure in Cl Tau. It is interesting to note that its /3 is close to 
0, which implies that the mass might be underestimated. On the other hand, as discussed 
in previous studies (e.g., Kwon et al. 2011), HL Tau has a gravitationally unstable region 
between 50 and 100 AU, where we marginally detected a possible substructure and the 
ALMA SV data present a wide gap. Also, DL Tau has a region around 100 AU with the Q 
values very close to the gravitational instability regime. ALMA observations toward these 
objects will provide a more clear view. 
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Fig. 1.— Protoplanetary disk images in A = 1.3 mm continuum. The color scales are in units 
of Jy beam“^, and the synthesized beams and the 100 AU bars are marked at the bottom of 
each panel. 
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Fig. 2.— uv amplitude plots of targets with their best htting models for the accretion disk 
model in solid lines and the power-law disk model in dashed lines. (5 values are calculated 
from the observational data simply assuming optically thin and Rayleigh-Jeans approxima¬ 
tions. The (5 results without these approximations are shown in the modeling. Open squares 
and triangles are for A = 1.3 mm and A = 2.7 mm data, respectively. 
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Fig. 3.— Posterior distributions of accretion disk model parameters fitting 6 disk data. 
Target names are indicated in the panels, the second column from the left. The M^isk is 
in units of M©, Rin and Rout are in AU, and inclination and position angles are in degree. 
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Fig. 4. 


Posterior distributions of power-law disk model parameters fitting 6 disk data. 
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Fig. 5.— Disk continuum, model, and residual maps at A = 1.3 mm and 2.7 mm for the ac¬ 
cretion disk models. From the left, individual target observational images, best models, and 
residual maps, and the upper and lower panels are at A = 1.3 mm and 2.7 mm, respectively. 
Note that the A = 1.3 mm images are the same as the color images of Figure [Tj The images 
are re-centered by the phase-center offsets listed in Table and the contour levels are ± 2.5, 
4.0, 6.3, 10, 16, 25, 40, and 100 times individual target a listed in Table The synthesized 
beams are marked in the bottom right corner of the observational images. 
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Fig. 6.— Relationship between properties: the mid-plane density gradient s vs. /9, M^isk vs. 
s, the flare index h vs. Rc, and Mdisk vs. f3 from the npper left clockwise. 
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Fig. 7.— Disk vertical structure parameter plots. The posterior is normalized by the peak 
and red indicates high posterior regions. 
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Fig. 8.— ALMA SV data of HL Tau in Band 7 (343.5 GHz) with an angular resolution 
of 0.03" X 0.02" (—176°) in color scales, overlaid with the residuals of our CARMA data at 
A = 1.3 mm after subtraction of the disk model in contours. The green and gray contours 
indicate positive and negative levels of 2, 3, and 4 times 0.75 mJy beam“^, respectively. The 
image in the right has an additional contour of 0.75 mJy beam“^ and the region out of the 
disk less than 3 mJy beam“^ has been masked. The coordinates are in epoch of 2014 Oct. 
31, when the ALMA data were taken, and our CARMA data were shifted by the proper 
motion. Note that the residuals nicely match the gaps and rings. 
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Fig. 9.— Comparisons of temperature distributions particularly for HL Tau in ALMA SV 
data and our modeling. Radial distributions of surface (T^, upper straight lines) and mid¬ 
plane (Tm, lower curves) temperatures are shown for each of the best htting accretion disk 
models. The black straight line indicates the power-law temperature distribution outlining 


the brightness temperature distribution of the ALMA SV data (Partnership et ah 2015): 
Tb[K] = 60(i?/20AU)“°'®^. The lower limit of dust temperatures in our modeling is 7 K. 
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Fig. 10.— The Toomre Q parameter along radius in the accretion disk model. 
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Table 1: Reference positions and stellar parameters 


Targets 

Positions (32000)'^ 
h m s ° " 

T.//” 

M 

[LO] 

[MO] 

EW(10/rm)® 

ni3-3i® 

Cl Tau 

04 33 52.000 

-b22 50 30.20 

4060 

1.40 

0.70 

2.55 

-0.17 

DL Tau 

04 33 39.076 

-b25 20 38.14 

4060 

1.12 

0.70 

0.51 

-0.77 

DO Tau 

04 38 28.583 

-b26 10 49.85 

3850 

2.70 

0.70 

0.94 

-0.13 

FT Tau 

04 23 39.178 

-b24 56 14.30 

3890^ 

0.38 

0.70 

1.86 

-0.46 

Haro 6-13 

04 32 15.410 

-b24 28 59.97 

3850^ 

2.11^ 

0.55 

3.87 

0.38 

HL Tau 

04 31 38.471 

-bl8 13 58.11 

4000^ 

8.30^ 

0.55 

- 

- 


^Observation phase center positions. Note that the actual disk center positions are offset as given in Table 


Kenyon & Hartmann (1995) 


'^Bolometric luminosities of Kenyon & Hartmann (1995) are adopted. 

'^These values are reasonably assumed, as they are not well known. Note that they are not a sensitive 
parameter for modeling. 


Furlan et al. (20091 


^It is assumed as same as or similar to Andrews & Williams (2005). 
^The values used in Kwon et al. (2011) have been employed. 






















Table 2: Protoplanetary disk observations 


Targets 

Obs. Dates 

Array 

Wavelengths 

Flux cal. 

Gain cal.. 

Cl Tau 

2007 Nov. 25 

B 

1 

MWC349 (1.9) 

0530-M35 (2.4), 3C84 (4.2) 

Cl Tau 

2007 Sep. 18 

C 

1 

Mars 

0530-M35 (3.0) 

Cl Tau 

2008 Jun. 13 

D 

1 

Uranus 

0530-M35 (1.6), 3C111 (1.7) 

Cl Tau 

2008 Jan. 18 

B 

3 

MWC349 (1.3) 

0530-M35 (3.8), 3C84 (7.1) 

Cl Tau 

2008 Oct. 29 

C 

3 

3C273 (14) 

0530-M35 (2.0) 

Cl Tau 

2009 Mar. 30 

D 

3 

Uranus 

0510-M80 (1.1), 3C111 (2.8) 

DL Tau 

2008 Dec. 11 

B 

1 

3C84 (5.0) 

0510-M80 (0.9) 

DL Tau 

2008 Dec. 14 

B 

1 

3C84 (5.0) 

0510-M80 (0.9) 

DL Tau 

2009 Mar. 11 

D 

1 

Uranus 

0510-M80 (0.9), 3C111 (2.3) 

DL Tau 

2009 Jan. 21 

A 

3 

3C454.3 (10) 

0510-M80 (1.3), 3C111 (3.8) 

DL Tau 

2009 Jan. 22 

A 

3 

3C84 (9.0) 

0510-M80 (1.3) 

DL Tau 

2008 Oct. 21 

C 

3 

Uranus 

0530-M35 (2.0), 3C111 (6.9) 

DL Tau 

2008 Oct. 29 

C 

3 

Uranus 

0530-M35 (2.0), 3C111 (6.9) 

DO Tau 

2007 Nov. 27 

B 

1 

MWC349 (1.9) 

0510-M80 (0.6) 

DO Tau 

2009 Nov. 4 

C 

1 

Uranus 

3C111 (1.2), 0510-M80 (0.5) 

DO Tau 

2009 Mar. 27 

D 

1 

Uranus 

0510-^180 (0.8), 3C111 (2.0) 

DO Tau 

2008 Oct. 18 

C 

3 

3C84 (9.0) 

0530-M35 (2.0) 

DO Tau 

2008 Oct. 28 

C 

3 

Uranus 

0530-M35 (2.0), 3C111 (6.9) 

DO Tau 

2008 Oct. 28 

C 

3 

3C84 (9.0) 

0530-M35 (2.0) 

FT Tau 

2011 Jan. 5 

B 

1 

Uranus 

3C111 (3.2), 0336-^323 (0.9) 

FT Tau 

2008 Oct. 15 

C 

1 

Uranus 

0357-^233 (0.35) 

FT Tau 

2009 Feb. 15 

A 

3 

Uranus 

0510-^180 (1.3), 3C111 (3.5) 

FT Tau 

2008 Oct. 26 

C 

3 

3C454.3 (20) 

0530-M35 (2.0), 3C111 (6.9) 

Haro 6-13 

2008 Jan. 12 

B 

1 

3C454.3 (11) 

0530-M35 (2.2) 

Haro 6-13 

2008 Apr. 25 

C 

1 

Uranus 

0530-^135 (2.0), 3C111 (1.9) 

Haro 6-13 

2010 Jan. 29 

A 

3 

Uranus 

3C111 (2.4), 0336-^323 (2.1) 

Haro 6-13 

2009 Jan. 4 

B 

3 

Uranus 

0510-^180 (1.3) 

Haro 6-13 

2009 Apr. 18 

C 

3 

Uranus 

0510-M80 (1.0), 3C111 (2.6) 

HL Tau 

2009 Jan. 17 

A 

1 


0510-M80 (0.9) 

HL Tau 

2009 Jan. 31 

A 

1 

3C84 

0510-^180 (0.9) 

HL Tau 

2009 Jan. 1 

B 

1 

3C84 (5.0) 

0510-^180 (0.9) 

HL Tau 

2007 Nov. 4 

C 

1 

Uranus 

0530-^135 (2.6) 

HL Tau 

2008 Feb. 16 

B 

3 

Uranus 

0530-^135 (4.0), 3C111 (6.7) 

HL Tau 

2008 Oct. 27 

C 

3 

Uranus 

0530-^135 (2.0) 



Table 3: Reduced data sets. The total fluxes are estimated as a flux sum over the target region beyond 2a levels and 
their uncertainties are statistical, which do not include the absolute flux calibration uncertainties. The ARA and ADec 
are offsets of target positions from the phase centers. 


Targets 

Arrays 

Freq. 

[GHz] 

uv coverage 
[kA] 

Beam 

["x" (PA°)] 

RMS 

[mJy beam“^] 

Total flux 
[mJy] 

Gaussian fit sizes 
["x" (PA°)] 

ARA 

n 

ADec 

n 

Cl Tan 

BCD 

229 

6.8—726.5 

0.45 X 0.31 

(-81) 

1.8 

97.8 ± 5.7 

0.82 X 0.64 (-14) 

0.24 

-0.20 


BCD 

113 

3.9—363.8 

0.95 X 0.69 

(89) 

0.55 

28.3 ± 1.5 

1.04 X 0.92 (13) 



DL Tan 

BD 

229 

8.5—723.8 

0.49 X 0.29 

(82) 

0.97 

194.2 ± 4.8 

1.02 X 0.78 (59) 

0.05 

-0.26 


AC 

113 

7.4—727.6 

0.31 X 0.27 

(-83) 

0.57 

35.0 ± 2.4 

1.14 X 1.11 (9) 



DO Tan 

BCD 

229 

7.0—616.6 

0.48 X 0.42 

(80) 

1.2 

125.4 ± 3.5 

0.59 X 0.47 (89) 

0.07 

-0.58 


C 

113 

8.3—143.8 

1.75 X 1.19 

(83) 

0.67 

30.8 ± 1.3 

point source 



FT Tan 

BC 

229 

17.0—620.0 

0.44 X 0.29 

(-71) 

1.1 

104.1 ± 3.3 

0.46 X 0.37 (-52) 

0.20 

-0.30 


AC 

113 

4.1—727.6 

0.31 X 0.28 

(-9) 

0.53 

28.8 ± 2.0 

1.27 X 0.85 (72) 



Haro 6-13 

BC 

229 

12.3—616.6 

0.54 X 0.39 

(-54) 

2.3 

111.2 ± 4.9 

0.34 X 0.28 (-7) 

0.15 

-0.59 


ABC 

113 

8.0—727.6 

0.34 X 0.29 

(-84) 

0.48 

28.7 ± 1.0 

0.26 X 0.20 (-69) 



HL Tan 

ABC 

229 

15.5—1452.0 

0.17 X 0.13 

(85) 

0.75 

685.0 ± 6.7 

0.85 X 0.64 (-46) 

-0.75 

-0.74 


BC 

112 

8.0—359.0 

1.01 X 0.67 

(73) 

0.93 

118.4 ± 2.6 

0.84 X 0.71 (-23) 





Table 4: Disk fitting results of free parameters. 


Viscous accretion disk model 

Targets 

s 

/3 

^disk 

[Mo] 

Rin 

[AU] 

Rc 

[AU] 

e^ 

n 

PA 

n 


Rt 

[AU] 

Cl Tau 

1.67±0.24 

0.049±0.027 

0.018±0.002 

4.2±2.6 

119.6±6.1 

61.0±2.9 

17.6±2.6 

(0.72) 

(57.4) 

DL Tau 

1.22±0.14 

0.408±0.022 

0.041±0.005 

1.3±1.1 

112.2±2.4 

-38.6±1.3 

55.2±2.0 

(0.27) 

(54.8) 

DO Tau 

1.38±0.21 

0.029±0.018 

0.014±0.001 

6.7±1.6 

54.2±1.1 

-32.5±2.0 

89.7±3.8 

(-0.17) 

(27.6) 

FT Tau 

1.85±0.24 

0.384±0.024 

0.021±0.002 

8.09±0.98 

55.2±2.0 

33.8±1.7 

136.1±2.9 

(0.60) 

(26.5) 

Haro 6-13 

2.16±0.21 

-0.152±0.024 

0.012±0.001 

1.6±1.3 

28.9±2.1 

41.5±3.3 

167.2±5.8 

(1.06) 

(14.8) 

HL Tau 

0.7483±0.0087 

0.6745±0.0069 

0.105±0.001 

8.78±0.19 

80.20±0.34 

41.19±0.22 

135.38±0.34 

(-0.20) 

(40.9) 

Power-law disk model 

Targets 

s 

/3 

^disk 

Rin 

Rout 

e^ 

PA 

p° 





[Mo] 

[AU] 

[AU] 

n 

n 



Cl Tau 

2.15±0.25 

0.036±0.027 

0.017±0.002 

6.5±2.0 

195±11 

62.5±2.6 

17.7±2.5 

(1.05) 


DL Tau 

2.08±0.23 

0.387±0.021 

0.029±0.003 

8.18±0.89 

152.6±2.3 

-38.4±1.3 

55.6±2.0 

(0.53) 


DO Tau 

1.98±0.23 

0.021±0.018 

0.013±0.001 

7.3±1.5 

75.4±2.0 

-32.8±2.0 

90.7±3.5 

(0.43) 


FT Tau 

2.80±0.15 

0.341±0.024 

0.015±0.002 

8.70±0.27 

98.2±5.1 

33.2±1.7 

136.0±3.0 

(1.25) 


Haro 6-13 

2.76±0.21 

-0.164±0.024 

0.011±0.002 

1.54±0.65 

61.9±5.2 

43.5±2.9 

167.3±5.9 

(1.81) 


HL Tau 

2.14±0.01 

0.615±0.006 

0.067±0.001 

8.70±0.02 

109.71±0.23 

40.9±0.20 

135.7±0.34 

(0.59) 



Note. — The presented uncertainties are statistical fit errors. As discussed in Table 1 of Kwon et al. 


(2011 , the systematic uncertainties are about 0.07 in s due to the flux calibration errors over different 


configuration data, 0.25 in j3 due to the flux calibration uncertainties over wavelengths, and 10% in Mdisk 
due to the flux calibration uncertainty at a dominant wavelength. The systematic uncertainties of Rim Rc 
and Rout depend on the angular resolution. The Rm uncertainty is also affected by the data sensitivity. 


“For computing this (s — h), the best h values of the accretion disk model are adopted: 0.95 (Cl Tau, DL 
Tau, HL Tau), 1.55 (DO Tau), 1.25 (FT Tau), and 1.10 (Haro 6-13). It is not intended to provide the best 
fitting result for surface density distributions in this paper. These are only for a rough comparison with 
other studies. 


^The best h values of the power-law disk model are adopted: 1.10 (Cl Tau), 1.55 (DL Tau, DO Tau, FT Tau, 
HL Tau), and 0.95 (Haro 6-13). 

















Table 5: Disk model comparison. 


Targets 

ln(P(D 1 Ha)) 

ln(P(D 1 Hp)) 

In(iF) 

Preferable model 

Cl Tau 

-1221129.9 

-1221128.1 

0.8 

Comparable 

DL Tau 

-5336360.9 

-5336365.6 

-4.7 

Power-law 

DO Tau 

-1794530.0 

-1794521.3 

8.7 

Accretion 

FT Tau 

-5444060.7 

-5444044.3 

16.4 

Accretion 

Haro 6-13 

-3053420.6 

-3053420.1 

0.5 

Comparable 

HL Tau 

-7863636.4 

-7863346.8 

289.6 

Accretion 
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Table 6: Correlation coefficients between properties. HL Tan is not included in the relation¬ 
ships with EW(10 /im) and ni 3 _ 3 i and the bold-font anti/correlations are plotted in Figure 

El 


Properties 

/5 

^disk 

Rin 

Rc 

h 

EW 

ni3-3i 

s 

-0.77 

-0.84 

-0.41 

-0.47 

0.21 

0.93 

0.74 

/? 

1.00 

0.84 

0.52 

0.36 

-0.36 

-0.70 

-0.95 

^disk 


1.00 

0.44 

0.26 

-0.48 

-0.65 

-0.86 

Rin 



1.00 

-0.16 

0.34 

-0.23 

-0.13 

Rc 




1.00 

-0.59 

-0.44 

-0.63 

h 





1.00 

-0.25 

0.19 

EW(10 pm) 

ni3-3i 

1.00 

0.83 

1.00 
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Table 7: Disk mass accretion properties. 


Targets 

Mdisk 

[xlO”"^ M© year“^] 

Mdisk^ / Mdisk 
[xlO® year] 

a(10 AU) 

a(100 AU) 

Cl Tau 

1.2 (0.49-2.6) 

1.5 

0.43 

0.18 

DL Tau 

1.2 (0.69-1.2) 

3.5 

0.30 

0.098 

DO Tau 

6.8 (1.5-6.8) 

0.2 

0.37 

0.017 

FT Tau 

0.51 (0.028-1.6) 

4.0 

0.029 

0.014 

Haro 6-13 

1.0^ 

1.2 

0.018 

0.013 

HL Tau 

10 (3.0-51) 

1.0 

0.41 

0.079 


“These values come from Robitaille et al. (2007). Particularly, we adopt the SED fitting values, which are 
available for most of our targets: the best fitting value and the ranges in parentheses. 

^The accretion disk model estimates in Tableware used. 


^Assumed as it is not available. 










